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Abstract 

In this paper we propose a Lagrangian method for solving Lane-Emden equation which is 
a nonlinear ordinary differential equation on semi-infinite interval. This approach is based 
on a Modified generalized Laguerre functions Lagrangian method. The method reduces the 
solution of this problem to the solution of a system of algebraic equations. We also present 
the comparison of this work with some well-known results and show that the present solution 
is acceptable. 
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1. Introduction 

Recently, spectral methods have been successfully applied in the approximation of dif- 
ferential boundary value problems defined in unbounded domains. We can apply different 
approaches using spectral methods to solve problems in semi-infinite domains. 

The first approach is using Laguerre polynomials/functions [3-0]. Guo [_J suggested a 
Laguerre- Galerkin method for the Burgers equation and Benjamin-Bona-Mahony (BBM) 
equation on a semi-infinite interval. It is shown that the Laguerre-Galerkin approximations 
are convergent on a semi-infinite interval with spectral accuracy. He in [_j introduced a new 
family of generalized Laguerre polynomials and investigated various orthogonal projections. 
Shen [_| proposed spectral methods using Laguerre functions and analyzed elliptic equa- 
tions on regular unbounded domains. In [3] is shown that spectral-Galerkin approximations 
based on Laguerre functions are stable and convergent with spectral accuracy in the Sobolev 
spaces. Siyyam [_J applied two numerical methods for solving initial value problem differen- 
tial equations using the Laguerre Tau method. Maday, et al. j_| proposed a Laguerre type 
spectral method for solving partial differential equations. 
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The second approach is reformulating the original problem in a semi-infinite domain to a 
singular problem in a bounded domain by variable transformation and then using the Jacobi 
polynomials to approximate the resulting singular problem |f|. 

The third approach replacing the semi-infinite domain with [0, K\ interval by choosing 
K, sufficiently large. This method is named domain truncation [7|. 

The fourth approach of spectral method is based on rational orthogonal functions. Boyd 
defined a new spectral basis, named rational Chebyshev functions on the semi-infinite 
interval, by mapping to the Chebyshev polynomials. Guo et al. [9J introduced a new set 
of rational Legendre functions which is mutually orthogonal in L 2 (0, +00). They applied 
a spectral scheme using the rational Legendre functions for solving the Korteweg-de Vries 



equation on the half line. Boyd et al. applied pseudospectral methods on a semi-infinite 



interval and compared rational Chebyshev, Laguerre, and mapped Fourier sine. 



The authors of |llMl3j applied spectral method to solve nonlinear ordinary differential 
equations on semi-infinite intervals. Their approach was based on a rational Tau method. 
They obtained the operational matrices of derivative and product of rational Chebyshev, 
Legendre functions, then applied these matrices together with Tau method to reduce the 
solution of these problems to the solution of a system of algebraic equations. The authors 



of [14J also applied pseudospectral method based on rational Legendre functions to solve 
Lane-Emden equations. 

This paper is arranged as follows: 
In section [2] we describe Lane-Emden equation. In section [3] we describe the formulation of 
generalized Laguerre polynomials and modified generalized Laguerre functions required for 
our subsequent development. In section H] we detailed description of Lagrangian interpolants 
construction and properties then obtained the operational matrices of derivative of modi- 
fied generalized Laguerre functions to applied these matrices together with the Lagrangian 
method to reduce the solution of this problem to the solution of the system of algebraic equa- 
tion. Section summarizes the application of this method for solving Lane-Emden equation 
and a comparison is made with existing methods in the literature. The results show pref- 
erence of this method in comparison with the others. The conclusions are described in the 
final section. 

2. Lane-Emden equation 



In the study of stellar structure |15[ an important mathematical model described by the 



second-order ordinary differential equation 

xy" + 2y' + xg(y) = 0, x > 0, (1) 
arises, where g(y) is some given function of y. Among the most popular form of g(y) is 

g{y) = y m , (2) 

where m is a constant, subject to the conditions 



y(0) = 1, y'(0) = 0. 
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(3) 



This equation is standard Lane-Emden equation. It was first proposed by Lane [16| and 
studied in more detail by Emden 

The Lane-Emden equation describes a variety of phenomena in theoretical physics and 
astrophysics, including aspects of stellar structure, the thermal history of a spherical cloud 



of gas, isothermal gas spheres, and thermionic currents 15 



2.1. Standard Lane-Emden equation 

This equation is one of the basic equations in the theory of stellar structure and has 
been the focus of many studies 18 -2(|. This equation describes the temperature variation 



of a spherical gas cloud under the mutual attraction of its molecules and subject to the 
laws of classical thermodynamics. The polytropic theory of stars essentially follows out of 
thermodynamic considerations, that deal with the issue of energy transport, through the 
transfer of material between different levels of the star. We simply begin with the Poisson 
equation and the condition for hydrostatic equilibrium: 



dP GM(r 
= - p 

dr r 



P^, (4) 



dM(r) , . 



dr 



Airpr , (5) 



where G is the gravitational constant, P is the pressure at radius r, M(r) is the mass of a 
star at a certain radius r, and p is the density, at a distance r from the center of a spherical 
star. Combination of these equations yields the following equation, which as should be 
noted, is an equivalent form of the Poisson Equation. 

1 d fr 2 dP\ A „ 

-— -— = -AnGp. (6) 



r 2 dr \ p dr / 

From these equations one can obtain the Lane-Emden equation through the simple sup- 
position that the pressure is related to the density, while remaining independent of the 
temperature. We already know that in the case of a degenerate electron gas that the pres- 
sure and density are p ~ Ps, assuming that such a relation exists for other states of the star 
we are led to consider a relation of the following form: 

P = Kp 1+ ™, (7) 

where K and m are constants, at this point it is important to note that m is the polytropic 
index which is related to the ratio of specific heats of the gas comprising the star. Based 
upon these assumptions we can insert this relation into our first equation for the hydrostatic 
equilibrium condition and from this rewrite equation to: 



AnG 



-\> 



1 d ( 



r 2 dr V dr 



3 



where the additional alteration to the expression for density has been inserted with A rep- 
resenting the central density of the star and y that of a related dimensionless quantity that 
are both related to p through the following relation 

P = AlT- (9) 

where m is a constant. 

Additionally, if place this result into the Poisson equation, we obtain a differential equa- 
tion for the mass, with a dependance upon the polytropic index m. Though the differential 
equation is seemingly difficult to solve, this problem can be partially alleviated by the in- 
troduction of an additional dimensionless variable x, given by the following: 

r = ax, (10) 



A m 

AttG 



Inserting these relations into our previous relations we obtain the famous form of the Lane- 
Emden equation, given below: 



x 2 dx 



x 



dx 



Taking these simple relations we will have the Lane-Emden equation with g(y) = y r 



y" + -y' + y r 

x 



0. 



x > 0. 



(12) 



(13) 



The physically interesting range of m is < m < 5. Numerical and perturbation approaches 
to solve equation Eq. (pQ) with g{y) = y m and boundary conditions (|3]) have been considered 
by various authors. It has been claimed in the literature that only for m = 0, 1 and 5 
the solutions of the Lane-Emden equation (also called the polytropic differential equations) 
could be given in closed form. 

In fact, for m — 5, only a 1-parameter family of solutions is presented. The so-called 



generalized Lane-Emden equation of the first kind have been looked at in Goenner [2l| and 
Havas 22 . 



2.2. Methods have been used to solve Lane-Emden equations 

Recently, many analytical methods have been used to solve Lane-Emden equations, the 
main difficulty arises in the singularity of the equation at x — 0. Currently, most techniques 
in use for handling the Lane-Emden-type problems are based on either series solutions or 
perturbation techniques. 

Bender et al. [l8|, proposed a new perturbation technique based on an artificial parameter 
6, the method is often called 5-method. 
Mandelzweig et al. 



231 ] used Quasilinearization approach to solve Lane-Emden equation. 



This method approximates the solution of a nonlinear differential equation by treating the 



nonlinear terms as a perturbation about the linear ones, and unlike perturbation theories is 
not based on the existence of some kind of a small parameter. He showed that the quasi- 
linearization method gives excellent results when applied to different nonlinear ordinary 
differential equations in physics, such as the Blasius, Duffing, Lane-Emden and Thomas- 
Fermi equations. 

Shawagfeh 24j applied a nonperturbative approximate analytical solution for the Lane- 
Emden equation using the Adomian decomposition method. His solution was in the form 
of a power series. He used Pade approximants method to accelerate the convergence of the 
power series. 

In 25(|, Wazwaz employed the Adomian decomposition method with an alternate framework 
designed to overcome the difficulty of the singular point. It was applied to the differential 
equations of Lane-Emden type. Further in 26| he used the modified decomposition method 
for solving analytical treatment of nonlinear differential equations such as Lane-Emden equa- 
tion. The modified method accelerates the rapid convergence of the series solution, dramati- 
cally reduces the size of work, and provides the solution by using few iterations only without 
any need to the so-called Adomian polynomials. 

Liao 27| provided a reliable, easy-to-use analytical algorithm for Lane-Emden type equa- 
tions. This algorithm logically contains the well-known Adomian decomposition method. 
Different from all other analytical techniques, this algorithm itself provides us with a con- 
venient way to adjust convergence regions even without Pade technique. 
He ' 
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employed Ritz's method to obtain an analytical solution of the problem. By the 
semi-inverse method, a variational principle is obtained for the Lane-Emden equation, which 
he gave much numerical convenience when applied to finite element methods or Ritz method. 



Parand et al. |lll-ll3l] presented some numerical techniques to solve higher ordinary differ- 



ential equations such as Lane-Emden. Their approach was based on a rational Chebyshev 
and rational Legendre tau method. They presented the derivative and product operational 
matrices of rational Chebyshev and rational Legendre functions. 

These matrices together with the tau method were utilized to reduce the solution of these 



physical problems to the solution of systems of algebraic equations. Parand et al. 14] also 
applied pseudospectral method based on rational Legendre functions to solve Lane-Emden 
equations. 

Ramos 29M32| solved Lane-Emden equation through different methods. In 29| he presented 
linearization methods for singular initial- value problems in second-order ordinary differential 
equations such as Lane-Emden. These methods result in linear constant-coefficients ordinary 
differential equations which can be integrated analytical, thus yielding piecewise analytical 
solutions and globally smooth solutions. Later, he (30| developed piecewise-adaptive decom- 
position methods for the solution of nonlinear ordinary differential equations. Piecewise- 
decomposition methods provide series solutions in intervals which are subject to continuity 
conditions at the end points of each interval, and their adaption is based on the use of either 
a fixed number of approximants and a variable step size, a variable number of approximants 
and a fixed step size or a variable number of approximants and a variable step size. 
In ' 
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series solutions of the Lane-Emden equation have been obtained by writing this 
equation as a Volterra integral equation and assuming that the nonlinearities are sufficiently 
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differentiable. These series solutions have been obtained by either working with the original 
differential equation or transforming it into an ordinary differential equation that does not 
contain first-order derivatives. It has been shown that these approaches provide exactly 
the same solutions as those based on Adomian's decomposition techniques that make use 
of either a different differential operator that overcomes the singularity at x = 0, or a new 
dependent variable, and Liao's homotopy analysis technique. Series solutions to the Lane- 
Emden equation have also been obtained by working directly on the original differential 
equation or transforming it into a simpler one. 

Yousefi 33] presented a numerical method for solving the Lane-Emden equations as singular 
initial value problems. Using integral operator and convert Lane-Emden equations to inte- 
gral equations and then applying Legendre wavelet approximations. He presented Legendre 
wavelet properties and then utilized these properties together with the Gaussian integration 
method to reduce the integral equations to the solution of algebraic equations. 

Chowdhury et al. presented a reliable algorithm based on the homotopy-perturbation 
method (HPM) to solve singular IVPs of time-independent equations. they obtained the ap- 
proximate and/or exact analytical solutions of the generalized Emden- Fowler type equa- 
tions. This method is a coupling of the perturbation method and the homotopy method. The 
HPM is a novel and effective method which can solve various nonlinear equations. The main 
feature of the HPM is that it deforms a difficult problem into a set of problems which are 
easier to solve. In this work, HPM yields solutions in convergent series forms with easily 
computable terms. 



Aslanov [35] introduced a further development in the Adomian decomposition method to 
overcome the difficulty at the singular point of non-homogeneous, linear and non-linear 
Lane-Emden-like equations; and constructed a recurrence relation for the components of 
the approximate solution and investigated the convergence conditions for the Emden- Fowler 
type of equations. He improved the previous results on the convergence radius of the series 
solution. Recently, Dehghan and Shakeri 36] first applied an exponential transformation to 
the Lane-Emden equation to overcome the difficulty of a singular point at x = and solved 



the resulting nonsingular problem by the variational iteration method. Yildirim et al. [37 
presented approximate exact solutions of a class of Lane-Emden type singular IVPs prob- 
lems, by the variational iteration method. The variational iteration method yields solutions 
in the forms of convergent series with easily calculable terms. Bataineh et al. [38] presented 
a reliable algorithm based on HAM to obtain the exact and/or approximate analytical solu- 
tions of the singular IVPs of the Emden-Fowler type. The HAM, first proposed by Liao in 



his Ph.D. dissertation [27], is a promising method for linear and non-linear problems. HAM 
contains an auxiliary parameter h which provides us with a simple way to adjust and control 
the convergence region and the rate of convergence of the series solution. 
Marzban et al. [39] used a method based upon hybrid function approximations. He used 
the properties of hybrid of block-pulse functions and Lagrange interpolating polynomials 
together with the operational integration matrix for solving nonlinear second-order, initial 
value problems and the Lane-Emden equation. 
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3. Properties of modified generalized Laguerre functions 

This section is devoted to the introduction of the basic notions and working tools con- 
cerning orthogonal modified generalized Laguerre functions. More specifically, we presented 
some properties of modified generalized Laguerre functions, concerning projection process. 

The Laguerre approximation has been widely used for numerical solutions of differential 
equations on infinite intervals. L°(x) (generalized Laguerre polynomial) is the nth eigen- 
function of the Sturm-Liouville problem 



x^L a n (x) + (a + 1 - x)^L a n (x) + nL a n (x) = 0, 
xG(0,oo), n = 0,l,2,.... , 
with the normalizing condition: 



n + a 
n 



|iS(«) = o, |l?(«) = -i, 

where a > — 1. 

The generalized Laguerre polynomials can be defined with the following recurrence formula: 
L%(x) = l, L%(x) = l + a-x, 

nL*{x) = (2n - 1 + a - x^^x) - (n + a - l)L%_ 2 (x), n>2 
these are orthogonal polynomials for the weight function w a = x a e~ x . The generalized 



Laguerre polynomials satisfy the following relation [41 



n-l 



d x L«{x) = -Y,L a k {x), (14) 

fc=0 

where n > 1 and a > — 1. 

We define Modified generalized Laguerre functions (which we denote MGL functions) T^(x) 
as follows: 

r^(ar) = exp(-x/(2L))I%(x/L), L > and a > -1. (15) 

This system is an orthogonal basis [l, 43] with weight function wl(x) = % and orthogonality 
property: 

^r a r a ^ - / £f> + a + l) \ 



L 2 n! 



where 5 nm is the Kronecker function. Boyd 0, 44] offered guidelines for optimizing the 
map parameter L where L > is the scaling parameter. On a semi-infinite domain, there is 
always a parameter that must be determined experimentally. 
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Numerical results deponed smoothly on constant parameter L, and therefore, are not very 
sensitive to L because the dError/dL = at the minimum itself, so the error varies very 
slowly with L around the minimum. A little trial and error is usually sufficient to find a 
value that is nearly optimum. In general, there is no way to avoid a small amount of trial and 
error in choosing L when solving problems on an unbounded domain. Experience and the 
asymptotic approximations of [44| can help, but some experimentation is always necessary 
as he explain in his book jjj. 



4. Lagrangian interpolation 

In this section we detailed description of Lagrangian interpolants construction and prop- 
erties, then we introduce the Lagrangian interpolation of modified generalized Laguerre 
functions and its operational matrices of derivative is obtained. 

Let {L%} <k<n is a basis in the space P n of polynomials of degree at most n. When n 
distinct points are given, another basis in P n is generated in a natural way. This is the basis 
of Lagrange polynomials with respect to the prescribed points. An element of the basis 
attains the value 1 at a certain point and vanishes in the remaining n — 1 points. 

Let us analyze first the generalize Laguerre polynomial case. We have the set of the 
Lagrange polynomials in P n relative to the n points rjk, < k < n — 1, i.e., the zeroes of L". 

The elements of the basis are denoted by < j < n — 1 (denoted Lagrangian 
interpolants). These polynomials in P n are uniquely defined by the conditions 

«s{m) = \\ t %= J- °<J< n ~ 1 - ( 16 ) 

[0 if i^j 

They actually form a basis because any polynomial p6P n can be written as follows: 

n-1 
3=0 

Therefore, p is a linear combination of the Lagrange polynomials. Such a combination is 
uniquely determined by the coefficients p{i]j), < j < n — 1. 
The following expression is easily proven: 

n-1 

TT x ~Vk 



II — 0<j'<n-l. 



'3 < X ~- 

For future applications, it is more convenient to consider the alternate expression 

if x ^ rjj, 




(17) 
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where < j < n — 1. Of course, we have 

L' a J 



lim 



x 



1, < J < 71 — 1. 



Lagrangian interpolants of generalized Laguerre polynomials (we denoted GLP) of order 



p at the Gauss-Radau-Laguerre quadrature points in M + is 41 



q(x) 



ig(0)' 



l,...,n, 



0. 



where r/j, j = 0, 1, 2, n are the n + 1 GLP-Radau points, 
derivative operator of GLP is: 

moreover for any polynomial p of degree at most n + 1, one gets: 



P 



'(^) = X^yP^i)- 



Funaro 0, 41] obtained derivative matrix of GLP(D n ): 



di 



ig(0) 



a+1 



i, j = 1, ...,n,i ^ j, 
i = j = l,...,n, 
i = 1, ...,n,j = 0, 
j = l, ... >r i, i = 0, 

i = j = 0. 



(19) 



The second derivative operator is obtained either by squaring D n either by evaluating 



^L°fa)((l-q+r ? ,)fa-^)-2r ?I ) 
fa -a) 2 _ n—1 

2(n+Q+l)Lg(0) 



n(n— 1) 
(a+l)(a+2) 



I, j = l,...,7l,l ^ j, 

i = l,...,n,j = 0, 
j = l,..., n ,2 = 0, 

i = j = 0. 



(20) 



Laguerre polynomials are not suitable for computations 41], and also Lagrangian interpo- 
lation of Laguerre polynomials is not suitable for solving some differential equations, such 
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as Lane-Emden equations because of their boundary conditions. So we use Lagrangian in- 
terpolation of MGL functions. At first we must find Lagrangian interpolants and derivative 
operators of MGL functions. 

Let T"(x) = e~ x / 2 L%(x) and by substitution of rji with x we have, 



X=Vi 



(21) 



X=T]i 



Lemma 1. Lagrangian interpolant ofT^(x) = e x ' 2 L c r 



\x) is, 



Q{x)=Q{x)- 



-x/2 



-Vi/2 ' 



(22) 



which £j(x) are Lagrangian interpolants of Laguerre polynomials. 

Proof: Suppose 7 xr "^ is Lagrangian interpolant of F^(x), using relation (T2TT) and 

x Vj 

L'Hopital's rule we can find constant 7, 

lim 7 ^1M = 7 lim (r-( x ) + x4-K{x)) = 1 

x-^rij x — T]j x->rij ax 

1 

^7 = 



so 



tAx) 



1 



='/J 



xe x / 2 L^{x) 



V^/ 2 £L-(x) 



X — T)j 



(23) 



and with comparison of ( Fl8|) and ( 123]) Lemma 1 is proved. 

By Eq. (T221 derivative operator of T°(x) = e~ x ' 2 L%(x) (we denote by _D n ) is: 



da da - 



-mil 



*3 ~ <*V e - m /2 



l/25 { 



i,j = 0,...,n. 



(24) 



Where matrix dij is defined in Eq. ( 119]) . As pointed out in before the second derivative 
operator is obtained either by squaring D n either by evaluating i^ n {r]i). For evaluating 
L n (rji) we can use the following relation: 



1 



- Vj /2 



4°v (l ' r ( 1, 2 



-m/2 



+ 



i, j = 0, 



n, 



(25) 



and iiPirjj) is defined in Eq. (TSUI) . It is obvious that MGL functions is T"(x/L), so La- 
grangian interpolant of MGL functions is ^(x/L), and derivative operators can be obtained 
easily. 
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4-1. Function Approximation 

We define the interpolant approximation of y(x) by 

n 
3=0 



(26) 



Where Pj(x) is defined in Eq. (123]) . The bj's are the expansion coefficients associated with 

the family {^(x/L)}. A semilogarithmic plot of \bj\ versus j is also useful to determine a 
good choice of L when the exact solution for y(x) is unknown. One can run the code for 
several different L and then plot the coefficients from each run on the same graph. The best 
L is the choice that gives the most rapid decrease of the coefficients 0. 
Therefore, 



I ny('3j)=b j , j = 0,...,n. 



(27) 



where Q^- = Lrjj are the zeroes of t n j{x/L). So derivative operator of MGL functions is (we 
denote by Dl u ): 



and the second derivative operator is 



L 2 



(2? 



(29) 



The relationship between the derivative 4-I n y(x) and I n y{x) at the collocation points Q\j 
0, ...,n can be obtained by differentiation. The result is as: 

, < i < n, 



h i di 

3=0 



(30) 



and 



d?x 



3=0 



0<i<n. 



(31) 



5. Solving Lane-Emden equation 

To apply Lagrangian interpolant of MGL functions to the Lane-Emden Equation intro- 
duced in Eq. ([[]) and Eq. (j2J) with boundary conditions Eq. (j3J), at first with Eq. f l2"B"j) we 
expand y(x), as follows: 

n 

I n y(x) = J2 b /j( x / L ), 

3=0 



11 



with a = 1 for Lagrangian interpolant of MGL functions. 

To find the unknown coefficients b/s, we substitute the truncated series into the Eq. (pQ) 
with g(y) introduced in Eq. ([2]) and boundary conditions in Eq. ([3]). So we have: 



^ E + 2 E h R + = °> 

J'=0 j=0 



l,...,n 



1 



(32) 



>o : 

n 



E M j 

3=0 



■oj 



0. 



(33) 



We have n—1 equations in Eq. ( 132)1 . that generates a set of n + 1 nonlinear equations with 
boundary equations in Eq. (133)1 . 

Table CD shows the approximations of y(x) for standard Lane-Emden with m = 3 obtained 
by the method proposed in this paper for n = 7 and L — 1, and those obtained by Horedt 

0- 

Table [2] shows the comparison of the first zero of y , between the present method, rational 
Legendre pseudospectral method (l4(, method in [39j, Pade approximation used by 18] and 
exact values reported in 45j for m = 2, 3, 4, respectively. 

Figure [T] shows the result graph of Lane-Emden for n = 6, m = 2, 3, 4. 



6. Conclusions 



The Lane-Emden equation describes a variety of phenomena in theoretical physics and 
astrophysics, including aspects of stellar structure, the thermal history of a spherical cloud 



of gas, isothermal gas spheres, and thermionic currents 15 



Since then the equation has been a center of attraction. The main problem in this direction 
is the accuracy and range of applicability of these approaches. The fundamental goal of 
this paper has been to construct an approximation to the solution of nonlinear Lane-Emden 
equation in a semi-infinite interval which has singularity at x = 0. A set of Laguerre 
functions are proposed to provide an effective but simple way to improve the convergence 
of the solution by Lagrangian method. Through the comparisons among the exact solutions 
of Horedt 



45 



the approximate solutions of Bender 18], recent good results in [14J and the 
current work, it has been shown that the present work has provided more accurate solutions 
for Lane-Emden equations, by taking to account N, it seems that the present approach 
gives the good result with small N. 
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Table 1: Approximation of y(x) for present method, solutions of Horedt |45| for to = 3 



X 


Present method 


0.000 


1.000000 


0.100 


0.998323 


0.500 


0.959821 


1.000 


0.855057 


5.000 


0.110820 


6.000 


0.043718 


6.800 


0.004165 


6.896 


0.000035 



solutions of Horedt [45 



1.000000 
0.998336 
0.959839 
0.855058 
0.110820 
0.043738 
0.004168 
0.000036 



Table 2: Comparison the first zero of y, between the present method, rational Legendre pseudospectral 



method [14J, method in [39(, Pade approximation used by [18| and exact values reported in [45| for to = 2, 3 
and 4 







Present 


method in [14] 


method 






m 


N 


method 


with iV = 75 


in [39] 


Bender 


Exact value 


2 


6 


4.352875 


4.35287108 


4.352875 


4.3603 


4.35287460 


3 


7 


6.896849 


6.89684862 


6.89685 


7.0521 


6.89684862 


4 


6 


14.971546 


14.9715463 


14.971546 


17.967 


14.9715463 
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Figure 1: Lane-Emden equation graph obtained by present method 
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